package org.platform.hrv.lib.algorithm;

import org.apache.commons.math3.complex.Complex;
import org.apache.commons.math3.stat.descriptive.moment.StandardDeviation;
import org.apache.commons.math3.transform.DftNormalization;
import org.apache.commons.math3.transform.FastFourierTransformer;
import org.apache.commons.math3.transform.FastHadamardTransformer;
import org.apache.commons.math3.transform.TransformType;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
import java.util.stream.Collectors;

public class HRVAnalysis {

    // 计算 SDNN（全部正常窦性心搏间期的标准差）
    public static double calculateSDNN(List<Double> rriList) {
        StandardDeviation stdDev = new StandardDeviation();
        return stdDev.evaluate(rriList.stream().mapToDouble(Double::doubleValue).toArray());
    }

    // 计算 RMSSD（相邻 RR 间期差值的均方根）
    public static double calculateRMSSD(List<Double> rriList) {
        double sumSquaredDiff = 0;
        int n = rriList.size();
        for (int i = 0; i < n - 1; i++) {
            double diff = rriList.get(i + 1) - rriList.get(i);
            sumSquaredDiff += diff * diff;
        }
        return Math.sqrt(sumSquaredDiff / (n - 1));
    }

    // 计算 NN50（相邻 RR 间期差值大于 50ms 的个数）
    public static int calculateNN50(List<Double> rriList) {
        int count = 0;
        int n = rriList.size();
        for (int i = 0; i < n - 1; i++) {
            double diff = Math.abs(rriList.get(i + 1) - rriList.get(i));
            if (diff > 50) {
                count++;
            }
        }
        return count;
    }

    // 计算 pNN50（NN50 占总心跳次数的百分比）
    public static double calculatepNN50(List<Double> rriList) {
        int n = rriList.size();
        int nn50 = calculateNN50(rriList);
        return (double) nn50 / n * 100;
    }

    // 计算频域指标 LF 和 HF
    public static void calculateFrequencyDomain(List<Double> rriList) {
        double[] rriArray = rriList.stream().mapToDouble(Double::doubleValue).toArray();
        FastFourierTransformer fft = new FastFourierTransformer(DftNormalization.STANDARD);
        Complex[] fftResult = fft.transform(rriArray, TransformType.FORWARD);

        double samplingFrequency = 1000.0 / averageRRInterval(rriList);  // 假设采样频率，根据平均 RR 间期计算
        double[] powerSpectrum = new double[fftResult.length];
        for (int i = 0; i < fftResult.length; i++) {
            powerSpectrum[i] = fftResult[i].abs() * fftResult[i].abs() / rriArray.length;
        }

        double lfPower = 0;
        double hfPower = 0;
        for (int i = 0; i < powerSpectrum.length; i++) {
            double freq = i * samplingFrequency / powerSpectrum.length;
            if (freq >= 0.04 && freq < 0.15) {
                lfPower += powerSpectrum[i];
            } else if (freq >= 0.15 && freq < 0.4) {
                hfPower += powerSpectrum[i];
            }
        }
        System.out.println("LF Power: " + lfPower);
        System.out.println("HF Power: " + hfPower);
    }

    // 计算平均 RR 间期
    public static double averageRRInterval(List<Double> rriList) {
        double sum = 0;
        for (Double rri : rriList) {
            sum += rri;
        }
        return sum / rriList.size();
    }

    public static void main(String[] args) {
//        List<Double> rriData = new ArrayList<>();
//        rriData.add(800.0);
//        rriData.add(820.0);
//        rriData.add(780.0);
//        rriData.add(850.0);
//        rriData.add(830.0);
//        rriData.add(790.0);
//        rriData.add(810.0);

        double[] doubles = {741,753,644,780,732,711,700,680,976,532,829,713,637,406,670,776,814,694,748,462,559,713,472,997,713,713,693,661,745,684,676,710,662,695,695,707,685,681,697,692,721,684,709,720,702,712,695,679,705,713,651,720,734,668,707,955,494,741,733,705,680,693,420,1171,685,528,796,533,1101,1140,845,582,713,527,1102,748,605,665,756,476,650,572,647,713,875,557,713,736,762,749,795,771,778,711,726,761,758,753,727,780,743,742,736,1143,736,686,714,718,745,748,742,787,794,792,730,768,750,755,719,705,740,754,738,718,747,747,750,737,742,740,784,733,703,738,720,731,680,699,718,735,753,745,740,708,710,744,722,724,722,740,745,761,723,717,741,722,705,674,727,707,733,700,730,729,673,770,737,755,745,758,710,731,766,486,670,673,525,673,688,729,710,729};
//        List<Double> rriData = Arrays.asList(
//                691.0, 794.0, 755.0, 601.0, 599.0, 1115.0, 661.0, 789.0, 1062.0, 440.0,
//                638.0, 642.0, 607.0, 1169.0, 578.0, 687.0, 679.0, 753.0, 715.0, 749.0,
//                639.0, 742.0, 789.0, 416.0, 789.0, 430.0, 694.0, 749.0, 829.0, 795.0,
//                833.0, 750.0, 830.0, 752.0, 871.0, 610.0, 781.0, 832.0, 789.0, 557.0,
//                814.0, 789.0, 789.0, 477.0, 1063.0, 618.0, 659.0, 691.0, 700.0, 639.0,
//                872.0, 811.0, 740.0, 588.0, 511.0, 666.0, 811.0, 811.0, 811.0, 895.0,
//                566.0, 811.0, 959.0, 637.0, 579.0, 793.0, 693.0, 811.0, 811.0, 671.0,
//                725.0, 799.0, 507.0, 858.0, 814.0, 811.0, 775.0, 771.0, 1077.0, 471.0,
//                465.0, 811.0, 811.0, 872.0, 861.0, 1076.0, 760.0, 702.0, 746.0, 663.0,
//                937.0, 811.0, 464.0, 1138.0, 811.0, 698.0, 621.0, 658.0, 581.0, 859.0,
//                670.0, 1047.0, 881.0, 811.0, 799.0, 811.0, 811.0, 723.0, 811.0, 623.0,
//                897.0, 1083.0, 811.0, 1085.0, 811.0, 811.0, 811.0, 811.0, 717.0, 782.0,
//                762.0, 455.0, 564.0, 1082.0, 781.0, 782.0, 782.0, 624.0, 775.0, 782.0,
//                782.0, 985.0, 567.0, 717.0, 798.0, 782.0
//        );
        List<Double> rriData = Arrays.stream(doubles).boxed().collect(Collectors.toList());

        // 计算时域指标
        double sdnn = calculateSDNN(rriData);
        double rmssd = calculateRMSSD(rriData);
        int nn50 = calculateNN50(rriData);
        double pnn50 = calculatepNN50(rriData);
        System.out.println("SDNN: " + sdnn);
        System.out.println("RMSSD: " + rmssd);
        System.out.println("NN50: " + nn50);
        System.out.println("pNN50: " + pnn50);

        // 计算频域指标
        calculateFrequencyDomain(rriData);
    }
}